\section{Iterative Löser für lineare Gleichungssysteme}
\subsection{Aufgabenstellung}
Gesucht ist $x\in \mathbb{R}^{n}$, so dass $Ax=b$ mit $A\in \mathbb{R}^{n\times n}, det(A)\not = 0,b\in \mathbb{R}^{n}$. Vorüberlegung:
\begin{enumerate}
	\item Liegt eine Bandmatrix mit Bandbreite $w$ vor, so ist der Aufwand zum exakten Lösen in $O(nw^{2})$
	\item Liegt eine dünnbesetzte Matrix vor, so ist der Aufwand einer Matrix-Vektor-Multiplikation in $O(n)$
\end{enumerate}

\begin{definition}
	[Definition III.1] Iterationsverfahren
	\begin{enumerate}
		\item Ein Lösungsverfahren zur Berechnung von $x$ mit $Ax=b$ heißt iterativ, falls ausgehend
		von einem Startwert $x^{0}$ eine Folge $x^{k}$ von Iterierten bestimmt wird.
		\item Ein Iterationsverfahren heißt konvergent, falls für alle $b$ und alle $x^{0}$ die Folge $x^{k}$ gegen $x$ konvergiert,
		d.h. $\lim_{k\rightarrow \infty} x^{k}=x$
		\item Ein Iterationsverfahren heißt konsistent, falls die Lösung $x$ Fixpunkt ist, d.h. falls $x^{k}=x\Rightarrow x^{k+1}=x$
	\end{enumerate}
\end{definition}

\subsection{Lineare Iterationsverfahren}

\begin{definition}
	[Definition III.2] Lineare Iterationsverfahren
	\\
	Ein Iterationsverfahren heißt linear, falls $x^{k+1}=Mx^{k}+Nb$ mit $M,N\in \mathbb{R}^{n\times n}$
\end{definition}

\begin{theorem} Konsistente lineare Iterationsverfahren
\\
Ein lineares und konsistentes Iterationsverfahren hat die Form:
$$x^{k+1}=x^{k}+N(b-Ax^{k})$$, d.h. $M=(Id-NA)$
\end{theorem}

\begin{theorem} Fehler
	\begin{enumerate}
		\item Für den Fehler $e^{k}=x-x^{k}$ gilt folgende Rekursion:
		$$e^{k+1}=M\cdot e^{k}=M^{k+1}\cdot e^{0}$$
		\item Ein lineares konsistentes Iterationsverfahren ist genau dann konvergent, falls $\rho(M)<1$ mit $\rho(M)$ der
		Spektralradius von $M$ (größter Eigenwert)
	\end{enumerate}
\end{theorem}
\begin{remark}
	$\rho(M)$  ist die Konvergenzrate
\end{remark}

\begin{definition}
	[Definition III.3] Lineare Iterationsverfahren
	\\
	Zerlege $A=L+D+U$
	\begin{description}
		\item[Richardson-Verfahren]
			$$x^{k+1}=x^{k}+w(b-Ax^{k}),\quad N=w\cdot Id \text{ mit $w$ Dämpfungsfaktor}$$
		\item[Jacobi-Verfahren]
			$$x^{k+1}=x^{k} + D^{-1}(b-Ax^{k}),\quad N=D^{-1},\quad M=-D^{-1}(L+U)$$
			$$x_{l}^{k+1}=x_{l}^{k}+\frac{(b_{l}-\sum_{i=1}^{n}a_{l,i}\cdot x_{i}^{k})}{a_{l,l}}$$
		\item[Gauß-Seidel-Verfahren]
			$$x^{k+1}=x^{k}+(D+L)^{-1}(b-Ax^{k}),\quad N=(D+L)^{-1},\quad M=-(D+L)^{-1}U$$
			$$x_{l}^{k+1}=x_{l}^{k}+\frac{(b_{l}-(\sum_{i=1}^{l-1}a_{l,i}\cdot x_{i}^{k+1} + \sum_{i=l}^{n}a_{l,i}\cdot x_{i}^{k}))}{a_{l,l}}$$
		\item[SOR-Verfahren]
			$$x^{k+1}=x^{k} + w(D+wL)^{-1}(b-Ax^{k}),\quad N=w(D+wL)^{-1}$$
			$$M=Id-w(D+wL)^{-1}A$$
	\end{description}
\end{definition}

\begin{remark}
	\begin{enumerate}
		\item Die Verfahren sind wohldefiniert, falls $a_{i,i}\not = 0$. Dies gilt z.B. für $A$ positiv definit oder $A$ diagonal dominant
		\item Das Richardson-Verfahren und das Jacobi-Verfahren sind unabhängig von der Indexanordnung. Das Gauß-Seidel- und das
		SOR-Verfahren sind abhängig von der Indexanordnung.
		\item Das Gauß-Seidel- und das SOR-Verfahren existieren auch in einer rückwärtigen bzw. symmetrischen Variante. Bei der 
		rückwärtigen Variante wird in der Definition von $N$ die Matrix $L$ durch $U$ ersetzt. Bei der symmetrischen Variante
		wird ein Standardhalbschritt und eine Rückwärtshalbschritt durchgeführt.
		\item Der Aufwand der Verfahren ist bei dünnbesetzten Matrizen $O(n)$ 
	\end{enumerate}
\end{remark}
\subsection{Konvergenzaussagen}
\begin{definition}
	[Definition III.4] Fehlerreduktion
	\\
	Man definiert:
	$$It(M)=-\frac{1}{ln(\rho(M))}$$
	als Maß für die Qualität des Verfahrens. Hierbei gibt $It(M)$ asymptotisch die Anzahl der Iterationsschritte für eine Fehlerreduktion
	von $\frac{1}{e}$ an.
\end{definition}

\begin{theorem}
	Konvergenzaussagen
	\\
	Sei die Matrix $A$ symmetrisch und positiv definit, dann konvergieren
	\begin{enumerate}
		\item Das Richardson-Verfahren für $w$ mit $0<w<\frac{2}{\lambda_{max}}$. Die optimale Konvergenzrate erhält man für
		$$w_{opt}=\frac{2}{\lambda_{max}+\lambda_{min}}\Rightarrow \rho(M_{opt})=1-\frac{2\lambda_{min}}{\lambda_{max}+\lambda_{min}}$$
		\item Das Jacobi-Verfahren konvergiert, falls $2D-A$ positiv definit ist.
		\item Das Gauß-Seidel-Verfahren konvergiert
		\item Das SOR-Verfahren koknvergiert für $0<w<2$
	\end{enumerate}
\end{theorem}

\begin{theorem}
	Konvergenzrate
	\\
	Sei $A$ symmetrisch positiv definit und $0<w<2$. Des Weiteren wird gefordert, dass die Eigenwerte der Matrix $B(z)=zD^{-1}L+\frac{1}{z}D^{-1}U$
	unabhängig von $z\in \mathbb{C}\setminus\{0\}$ sind. Dann gilt:
	\begin{enumerate}
		\item $\rho(M_{Jac})=:\beta<1$
		\item $\rho(M_{GS})=\rho(M_{Jac})^{2}$ (d.h. das Gauß-Seidel-Verfahren benötigt halb so viele Schritte wie das Jacobi-Verfahren)
		\item $\rho(M_{SOR,w})=
		\begin{cases}
		w-1 & w_{opt}\leq w<2\\
		1-w+\frac{1}{2}w^{2}\beta^{2}+w\beta\sqrt{1-w+\frac{w^{2}\beta^{2}}{4}} & 0 < w < w_{opt}
		\end{cases}$
	\end{enumerate}
	$w_{opt}=\frac{2}{1+\sqrt{1-\beta^{2}}}\geq 1$
\end{theorem}

\begin{remark}
	In der Regel ist $w_{opt}$ nicht berechenbar, da $\beta$ unbekannt ist. Ausweg: Adaptive Bestimmung einer Näherung $w_{opt}$.
\end{remark}

\subsection{Das konjugierte Gradientenverfahren (CG)}
Für symmetrisch positiv definite Matrizen $A$ gibt es ein nichtlineares Verfahren, welches im Vergleich zum
\begin{itemize}
	\item Gauß-Seidel deutlich schneller ist
	\item SOR mit $w_{opt}$ wird keine Zusatzinformation benötigt
\end{itemize}
Die Vorstufen des CG-Verfahrens sind
\begin{itemize}
	\item konjugierte Richtungsverfahren
	\item Gradientenverfahren
\end{itemize}

\begin{theorem}
	Sei $A$ symmetrisch positiv definit, so gilt: $Ax=b \Leftrightarrow F(x)\leq F(y),\ y\in \mathbb{R}^{n}$ 
	mit $F(y)=\frac{1}{2}y^{T}Ay-b^{T}y$
\end{theorem}

\subsubsection{Verfahren der konjugierten Richtungen}
\begin{definition}
	[Definition III.5] A-orthogonal
	\\
	Die Suchrichtungen $p^{0},p^{1},\ldots$ heißen A-orthogonal $\Leftrightarrow \langle p^{l},Ap^{k}\rangle=0$ mit $l\not = k$
\end{definition}

\begin{algorithm}
	\caption{Konjugierte Richtungen}
	\begin{algorithmic}
		\STATE  Seien $A$ symmetrisch positiv definit, $x^{0}$, und $p^{0},p^{1},\ldots,p^{n-1}$ A-othogonal
		\STATE $r^{0}=b-Ax^{0}$
		\FOR{$k=0$ to $n-1$}
			\STATE $\alpha_{k} = \frac{\langle p^{k},r^{k}\rangle}{\langle p^{k},Ap^{k}\rangle}$
			\STATE $x^{k+1} = x^{k} + \alpha_{k}p^{k}$
			\STATE $r^{k+1} = r^{k} - \alpha_{k}Ap^{k}$
		\ENDFOR
	\end{algorithmic}
\end{algorithm}

\begin{theorem}
	$x^{n}=x$
\end{theorem}
\begin{remark}
	Nachteile
	\begin{enumerate}
		\item Im ungünstigsten Fall passiert $(n-1)$-Mal nichts
		\item A-orthogonale Suchrichtungen sind teuer zu konstruieren
	\end{enumerate}
\end{remark}

\subsubsection{Gradientenverfahren}

\begin{algorithm}
	\caption{Gradientenverfahren}
	\begin{algorithmic}
		\STATE  Seien $A$ symmetrisch positiv definit, $x^{0}$
		\STATE $r^{0}=b-Ax^{0}, p^{0}=r^{0}$
		\FOR{$k=0$ to $\infty$}
			\STATE $\alpha_{k} = \frac{\langle p^{k},r^{k}\rangle}{\langle p^{k},Ap^{k}\rangle}$
			\STATE $x^{k+1} = x^{k} + \alpha_{k}p^{k}$
			\STATE $r^{k+1} = r^{k} - \alpha_{k}Ap^{k}$
			\STATE $p^{k+1}=r^{k+1}$
		\ENDFOR
	\end{algorithmic}
\end{algorithm}

\begin{remark}
	Beobachtung $x^{k+1}=x^{k}+\alpha_{k}(b-Ax^{k})$. Gleiche Struktur wie das Richardson-Verfahren. Aber $\alpha_{k}$ ist von $k$
	abhängig und damit ist es kein lineares Verfahren.
\end{remark}

\begin{theorem}
	Konvergenz
	\\
	$$||x^{k}-x||_{A}\leq \left( \frac{\kappa - 1}{\kappa +1}\right)^{k}\cdot || x^{0}-x||_{A}$$
	mit $\kappa = \frac{\lambda_{max}}{\lambda_{min}}$ und $||y||^{2}_{A}=\langle y, Ay \rangle$ ist die Energienorm.
\end{theorem}

\begin{description}
	\item[Vorteil] Das Gradientenverfahren konvergiert genauso gut wie das optimal gedämpfte Richardson-Verfahren ohne dessen Zusatzinformationen
					zu benötigen.
	\item[Nachteil] Minimierungseigenschaft bzgl. vorausgegangener Suchrichtungen geht i.d.R. verloren.
\end{description}

\subsubsection{CG-Verfahren}
Grundidee:
\begin{itemize}
	\item Baue im Lauf der Iteration ein System von $\dim$ unabhängigen A-orthogonale Suchrichtungen auf.
	\item Korrigiere die Suchrichtungen des Gradientenverfahrens
\end{itemize}

\begin{algorithm}
	\caption{CG}
	\begin{algorithmic}
		\STATE $A$ symmetrisch positiv definit, $x^{0}$, $b$
		\STATE $p^{0} = r^{0} = b- Ax^{0}$
		\FOR{$k=0$ to $\infty$}
			\STATE $\alpha_{k}=\frac{\langle r^{k},p^{k}\rangle}{\langle A p^{k},p^{k}\rangle}$
			\STATE $x^{k+1} =  x^{k} + \alpha_{k} p^{k}$
			\STATE $r^{k+1} = r^{k} - \alpha_{k}Ap^{k}$
			\STATE $\beta_{k} = - \frac{\langle r^{k+1},Ap^{k}\rangle}{\langle p^{k},Ap^{k}\rangle}$
			\STATE $p^{k+1} = r^{k+1} + \beta_{k}p^{k}$
		\ENDFOR
	\end{algorithmic}
\end{algorithm}

\begin{remark}
	Orthogonalitäten:
	\begin{itemize}
		\item $\langle r^{k+1}, p^{k}\rangle=0$
		\item $\langle p^{k+1},Ap^{k}\rangle=0$
	\end{itemize}
\end{remark}

\begin{theorem}
	Äquivalente Darstellung
	\\
	\begin{itemize}
		\item $\alpha_{k}=\frac{||r^{k}||_{2}^{2}}{||p^{k}||_{A}^{2}}$
		\item $\beta_{k} = \frac{||r^{k+1}||^{2}_{2}}{||r^{k}||^{2}_{2}}$
	\end{itemize}
\end{theorem}


\begin{theorem}
	Eigenschaften der Suchrichtungen
	\\
	Sei $x^{l}\not = x$ für $l=0,\ldots,k$ so gilt
	\begin{enumerate}
		\item $\langle r^{l},p^{j}\rangle=0\quad 0\leq j<l\leq k$
		\item $\langle r^{l},r^{j}\rangle=0\quad 0\leq j<l\leq k$
		\item $\langle p^{l},Ap^{j}\rangle =0 \quad 0\leq j<l\leq k$
	\end{enumerate}
\end{theorem}

\begin{definition}
[Definition III.6] Krylovraum $K_{l}$
\\
	$K_{l}:=K_{l}(r^{0},A):=span\{r^{0},Ar^{0},A^{2}r^{0},\ldots,A^{l-1}r^{0}\}$
\end{definition}

\begin{theorem}
	Suchrichtungen
	\\
	Für den Raum der Suchrichtungen $V_{k} :=span\{p^{0},p^{1},\ldots,p^{k-1}\}$ gilt: $V_{k}=K_{k}$ falls $x^{l}\not = x,l=0,\ldots,k$
\end{theorem}

\begin{theorem}
	Konvergenzaussage für das CG
	\\
	Sei $A$ symmetrisch positiv definit, so gilt für den Fehler $e^{k}$
	$$||e^{k}||_{A}\le \frac{2c^{k}}{1+c^{2k}}\cdot ||e^{0}||_{A}\quad c:=\frac{\sqrt{\kappa(A)}-1}{\sqrt{\kappa(A)}+1}$$
\end{theorem}

\begin{remark}
	Das CG-Verfahren konvergiert deutlich schneller als das Gradientenverfahren. Aber beide Verfahren konvergieren langsam, falls $\kappa(A)>>1$
\end{remark}

\subsubsection{Vorkonditionierte CG-Verfahren}
Sei $A$ symmetrisch positiv definit mit $\kappa(A)>>1$ und $B$ sei symmetrisch positiv definit mit $\kappa(AB)<<\kappa(A)$. Dann
ist
\begin{enumerate}
	\item $Ax=b \Leftrightarrow B^{\frac{1}{2}}AB^{\frac{1}{2}}y=B^{\frac{1}{2}}b\quad x=B^{\frac{1}{2}}y$
	\item Die Konvergenz des CG-Verfahren angewendet auf $B^{\frac{1}{2}}AB^{\frac{1}{2}}y=B^{\frac{1}{2}}b$ ist deutlich schneller als
	das CG-Verfahren angewandt auf $Ax=b$
\end{enumerate}

\begin{algorithm}
	\caption{Vorkonditioniertes CG-Verfahren}
	\begin{algorithmic}
		\STATE $A,B s.p.d.,b,x^{0}$
		\STATE $r^{0}=b-Ax^{0}$
		\STATE $p^{0}=Br^{0}$
		\STATE $\alpha_{0}=\langle r^{0},p^{0}\rangle$
		\FOR{$k=0$ to $size(A,1)$ if $\alpha_{k}\not = 0$}
			\STATE $v=Ap^{k}$
			\STATE $\lambda = \frac{\alpha_{k}}{\langle p^{k},v\rangle}$
			\STATE $x^{k+1}=x^{k}+\lambda p^{k}$
			\STATE $r^{k+1}=r^{k} - \lambda v$
			\STATE $z=Br^{k+1}$
			\STATE $\alpha_{k+1}=\langle r^{k+1},z\rangle$
			\STATE $p^{k+1} = z + \frac{\alpha_{k+1}}{\alpha_{k}}p^{k}$
		\ENDFOR
	\end{algorithmic}
\end{algorithm}

\begin{remark}
	Ein guter Vorkonditionierer erfüllt:
	\begin{enumerate}
		\item $\kappa(AB)=O(1)$
		\item $By\rightarrow $ Aufwand $O(n)$
	\end{enumerate}
\end{remark}